\documentclass[10pt,a4paper]{article}
\usepackage{graphicx}
\usepackage[english]{babel}
\usepackage{fullpage}
\usepackage{latexsym}
\usepackage{amssymb}
\usepackage[T1]{fontenc}
\usepackage[sc]{mathpazo}
\linespread{1.05}


\title{Applied Statistics \\ Presentation}
\author{Mieke Hiltermann and Lotte van den Berg \\ Utrecht University}

\newcommand{\mco}{\mathcal{O}}
\newcommand{\mcg}{\mathcal{G}}
\newcommand{\mcv}{\mathcal{V}}
\newcommand{\mce}{\mathcal{E}}
\newcommand{\mcn}{\mathcal{N}}
\newcommand{\real}{\mathbb{R}}



\begin{document}


\maketitle 
%\newpage

%\tableofcontents

\section{Introductie}

Paper \textit{`An SPC case study on stabilizing syringe lengths'} van L.A. Franklin en S.N. Mukherjee. 
Over de lengte van bepaalde spuitjes die gemaakt worden in een farmaceutisch bedrijf.\\

\noindent Tijdens het productieproces worden injectie spuiten gevuld met een enkele dosis van een bepaald medicijn. Eerst wordt het medicijn in glazen cartridges gedaan en geseald met een rubberen dop. In de fase die wordt bekeken met behulp van SPC worden the glazen cartridges in de plastic spuiten gedaan vervolgens wordt een stop op de spuitjes `geschoten'. Deze stop moet op een bepaalde lente van de spuit terecht komen. Deze lengte moet tussen de $4.290$ en de $4.980$ inch liggen want als de stop op minder dan $4.290$ inch ligt kan er teveel druk in de spuit ontstaat en dan kan (een deel van) het medicijn er al uit lopen; terwijl bij een te grote lengte de stop juist weer uit de spuit kan schieten. \\

\noindent Het artikel gaat over de implementatie van SPC-procedures om de de bovenstaande `tackingoperatie' te optimaliseren. 
De lengte van een spuitje de `spuitjeslengte' is de \textit{critical variable} en die gaan wordt bekeken met behulp van $\overline X$ en $R$ charts. \\

\noindent Opbouw presentatie:
\begin{itemize}
\item{capability study}
\item{$\overline X$ en $R$ charts}
\item{Monitoring en Adjustment stages}
\item{Statistical Quality Specialist}
\item{Wat was er aan de hand?}
\item{Conclusie}
\end{itemize}

\section{Capability Study}
Om de capability van het proces te testen werd de technicus gevraagd om de machine zo goed mogelijk in te stellen (!!! op gevoel).
Er werden 35 steekproeven genomen en er werd een capability-studie gedaan. \\

\noindent
Wat vooraf moet gaan aan de capabilitytest is een algemene controle van de data.
De gebruikte capabilitytest is gebaseerd op de aanname dat de data normaal verdeeld is. Dit moet eerst gecheckt worden.

\begin{verbatim}
paper <- read.table("paperfranklindata1.txt")

#informal check for outliers, Box and Whisker Plot:
boxplot(paper, horizontal = TRUE, main="Box-and-Whisker plot of initial data")

#Kernel density plot:
plot(density(paper$V1),main="Kernel density estimate of initial data",col="red",lwd=3)

#Normal probability plot:
qqnorm(paper$V1,main="Normal probability plot of initial data", pch=19,cex=1,fg="red")
qqline(paper$V1,lwd=3,col="blue",lty="dashed")
\end{verbatim}

\noindent
--> laten zien dmv boxplot; kernel density plot en normal probability plot \\
Als we de normality in R echter testen dmv een boxplot, kernel density plot en normal probability plot dan is het antwoord een beetje twijfelachtig:\\
- Er zijn geen outliers (alle data tussen $4.920$ en $4.980$)\\
- Linker kant van de densityplot is niet echt mooi, rechter kant wel.\\
- normal probabilityplot laat ook zien dat de bovenkwantielen er best mooi normaal uit zien, maar de onderkwantielen niet.\\
\\
Toch word besloten de capabilitytest te gebruiken: 

\begin{verbatim}{
#Confidence levels:
samplenumber <- rep(1:35, each=1)
paperdata <- qcc.groups(paper$V1, samplenumber)
paperqcc <- qcc(paperdata, type="xbar",spec.limits=c(4.920, 4.980))# 
process.capability(paperqcc, spec.limits=c(4.920, 4.980))}

\end{verbatim}
We zien dat het gemiddelde best dichtbij de gewenste waarde $4.950$ inch zit. Maar in de Boxplot zien we juist dat de driekwart van de data juist boven $4.950$ zit en de mediaan ongeveer bij $4.956$\\
\\
Er wordt dus de aanname gemaakt dat de data normaal verdeeld zijn terwijl je daar behoorlijk je twijfels bij kunt hebben. \\
\\
Daarnaast moet je natuurlijk wel om hieruit de control limits te halen enigszins ervan overtuigd zijn dat deze 35 datapunten in control waren.\\
\\
De capabilitytest laat zien dat het gemiddelde van de data op $\overline X=4.954$ zit, dit zit dichtbij de gewenste waar de van $4.950$. De waarden voor de capabilityindices zijn $C_p=1.2$ en $C_{pk}=1.02$. \\
\begin{verbatim}
> shapiro.test(paper$V1)

        Shapiro-Wilk normality test

data:  paper$V1 
W = 0.9393, p-value = 0.05346
\end{verbatim}
Lage p-value, maar net hoog genoeg om niet te verwerpen.\\
%\noindent UIt de plots hebben we gezien dat het process niet gecentreerd is, dus de waarde van $C_p$ zegt niet erg veel en we kunnen dus beter kijken naar de waarde van $C_{pk}$. 
Hieruit wordt geconcludeerd dat het proces ``minimaal capabel'' is en dat de goede lengtes coor de spuitjes gerealiseerd kan worden.\\
--> De minimale waardes voor $C_p$ en $C_{pk}$ waren eerst echter $1.33$ of $1.67$ en tegenwoordig worden zelfs vaak waarden van boven de $2.00$ `ge\"eist'. Als je hiernaar kijkt lijken deze waarden dus echt veel te laag. De vraag is dus of dit process wel capabel is.\\
\\
\noindent
De grafieken waarop deze data gebaseerd is zijn te zien in Figuur \ref{fig:PCA}.\\

\begin{figure*}[hbtp]  
\[
\begin{array}{ll}
\includegraphics[width=7cm]{paper-boxwhisker.jpg}&\includegraphics[width=7cm]{paper-densityplot.jpg}\\
\includegraphics[width=7cm]{paper-normalprobplot.jpg}&\includegraphics[width=7cm]{paper-pca.jpg}\\
\end{array}
\]
\caption{Process capability analysis}\label{fig:PCA}
\end{figure*}

% onderstaande wel of niet toevoegen?? 
% -----------------------------

\noindent We kijken nog even naar de proportie van non-conforming items. Omdat het proces niet gecentreerd is kunnen we dat als volgt berekenen.
$\frac{1}{2} (LSL + USL) = \frac{4.9200 + 4.9800}{2} = 4.9500 \leq \mu=4.9544 \leq 4.9800 = USL$. Dus de verwachte proportie non-conforming items is gelijk aan:
$$ 1- P(LSL < X < USL) = \Phi(-3(2C_p - C_{pk})) + \Phi(-3C_{pk})$$
\begin{verbatim}
> pnorm(-3*(2*1.195 - 1.02)) + pnorm(-3*1.02)
[1] 0.001126468
\end{verbatim}

% ---------------------------
\newpage
\section{$\overline X$ en $R$ charts}
Om de uiteindelijke gewenste control charts te krijgen werden er elke 15 minuten samples genomen met een grootte van 5, deze data komt uit het bestand \textit{`paperfranklindataallekollommen.txt'} (de eerste 15 rijen).\\
Deze data wordt gebruikt om de Upper en Lower Control Limits te berekenen. Dus eigenlijk wordt deze data gebruikt om te kijken of het process (phase II)in-control is, maar dan moet dat eerst waar zijn voor phase I. Daarom hebben we besloten hiervoor eers een Process Capability Indices uit te voeren., deze is te zien in figure \ref{fig:PCAphase1}.


\begin{figure*}[hbtp]  
\[
\begin{array}{ll}
\includegraphics[width=7cm]{paper-boxwhiskerphase1.jpg}&\includegraphics[width=7cm]{paper-densityplotphase1.jpg}\\
\includegraphics[width=7cm]{paper-normprobplotphase1.jpg}&\includegraphics[width=7cm]{paper-pcaphase1}\\
\end{array}
\]
\caption{Process capability analysis of phase I data}\label{fig:PCAphase1}
\end{figure*}

Het is erg twijvelachtig of deze data in-control is. Maar dit wordt wel aangenomen.\\
\\
\begin{verbatim}
> shapiro.test(paper2[1:15,7])

        Shapiro-Wilk normality test

data:  paper2[1:15, 7] 
W = 0.9579, p-value = 0.6553
\end{verbatim}

Om de upper en de lower control limits te bepalen hebben we eerst het gemiddelde nodig van de $R$ en de $\overline X$ charts.

% hier doen we dus geen capability studie voor .. nog wel doen !?! (inclusief normalitytesting)

\begin{verbatim}
> stats.R(paper2[1:15,2:6])
$statistics
    1     2     3     4     5     6     7     8     9    10    11    12    13 
0.014 0.031 0.042 0.042 0.014 0.018 0.016 0.020 0.056 0.036 0.016 0.032 0.012 
   14    15 
0.016 0.020 

$center
[1] 0.02566667

> stats.xbar(paper2[1:15,2:6])
$statistics
    1     2     3     4     5     6     7     8     9    10    11    12    13 
4.954 4.942 4.951 4.961 4.957 4.959 4.952 4.959 4.954 4.954 4.967 4.959 4.970 
   14    15 
4.963 4.970 

$center
[1] 4.958133
\end{verbatim}

\noindent Dus $E(X_i)=\mu=4.958133$ en $E(R_i)=0.02566667 $ voor $i=1,\ldots,15 $, ook hebben we de standaard deviaties nodig.
\begin{verbatim}
> sd.xbar(paper2[1:15,2:6])
[1] 0.01162417
> sd.R(paper2[1:15,2:6])
[1] 0.01162417
\end{verbatim}
dit resulteert in $\sigma=0.0116$ voor deze data. \\
\\
\noindent
Normaal gesproken worden de upper en lower control limits voor een $\overline X$ chart bepaald op de volgende manier:
\begin{eqnarray*}
(LCL,UCL)_{\overline X,\, normaal}&=&(\mu-\frac{3\sigma}{\sqrt n},\, \mu+\frac{3\sigma}{\sqrt n})\\
&=&(4.942538,\, 4.973728).\\
\end{eqnarray*}
In het artikel werd besloten de volgende waarden te pakken voor de control limits.
\begin{eqnarray*}
(LCL,UCL)_{\overline X,\, arikel}&=&(\mu-A_2(n) E(R_i),\, \mu+A_2(n) E(R_i))\\
&=&(4.958133-0.577*0.02566667,\, 4.958133+0.577*0.02566667)\\
&=&(4.9433, 4.9729)\\
\end{eqnarray*}
hierin is $A_2$ een Control Chart constante, deze is afhankelijk van de grootte van de rationele deelgroepen (in dit geval 5). 
Het is duidelijk dat het verschil tussen de verschillende control limits niet echt noemenswaardig is. Daarom gaan wij verder met de controllimits gegeven in het artikel.\\
\\\noindent
Normaal gesproken worden de control limits van de $R$ chart gezet op:
\begin{eqnarray*}
LCL_{R,\, normaal}&=&D_{0.001}(n)\frac{E(R_i)}{d_2(n)}=0.199*\frac{E(R_i)}{2.326}=0.004049728\\
UCL_{R,\, normaal}&=&D_{0.999}(n)\frac{E(R_i)}{d_2(n)}=5.484*\frac{E(R_i)}{2.326}= 0.0605142.\\
\end{eqnarray*}
In het artikel besloten ze echter om de volgende limieten te nemen. 
\begin{eqnarray*}
(LCL, UCL)_{R, \,artikel}&=&(D_3 E(R_i), D_4 E(R_i))\\
&=&(0.000*0.02566667, 2.115* 0.02566667)\\
&=&(0.000, 0.0543)\\
\end{eqnarray*}
Dit verschil is weer niet zo groot en we houden ons weer aan de limieten in het artikel.
%misschien een leuke toevoeging om ook R en X charts te maken met de 'normale'waarden...

Van deze data (die eerste 15 kolommen) en de limieten die we net hebben vastgesteld, kunnen we een $\overline X$ chart en een $R$ chart maken, deze zijn te zien in Figuur \ref{fig:control}\\
\begin{figure*}[hbtp]
\includegraphics[width=7cm]{paper-xbar.jpg} \includegraphics[width=7cm]{paper-R.jpg}\\
\caption{Control charts}\label{fig:control}
\end{figure*}
\begin{verbatim}
paper2 <- read.table("paperfranklindataallekollommen.txt")
paperqccX1 <- qcc(paper2[1:15,2:6], type="xbar", center=mu, limits=c(LCLXa,UCLXa))
paperqccR1 <- qcc(paper2[1:15,2:6], type="R", center=ER, limits=c(LCLRa,UCLRa))
\end{verbatim}
In deze charts is het zichtbaar dat het process nu al out-of-control is. Het gemiddelde van het proces ligt niet goed (het gemiddelde volgt een stijgende lijn) en de variatie is erg groot. Er wordt echter alleen gekeken of de data binnen de limieten blijft, dus hier wordt niks mee gedaan.

% aantal non-conforming items nog doen? 
% (bij gemiddelde en stijgende lijn refereren naar opmerking in conclusie van artikel: "8 in a row" rule ??)--> vind ik niet nodig

\section{The Process Monitoring and Adjustment Stages}

Hierna werden nog steeds elke 15 minuten metingen gedaan in groepen van 5. Deze metingen werden pas weer in een chart gezet nadat de 15 metingen gedaan waren. Dit zijn dus de metingen 16 tot en met 30. 
De upper en lower control limits die gebruikt worden voor deze data, zijn de limits van de eerste 15 metingen. Het resultaat is te zien in Figuur \ref{fig:control2}.
In dit figuur zijn ook oranje punten zichtbaar, deze komen als er meer dan 6 opeenvolgende punten boven of onder de center line komen. Dit zijn dus een soort van waarschuwingspunten.
\begin{figure*}[hbtp]
\includegraphics[width=7cm]{paper-xbar2.jpg} \includegraphics[width=7cm]{paper-R2.jpg}\\
\caption{Control charts}\label{fig:control2}
\end{figure*}
\begin{verbatim}{
paperqccX2 <- qcc(paper2[16:32,2:6], type="xbar", center=mu, limits=c(LCLXa,UCLXa))
paperqccR2 <- qcc(paper2[16:32,2:6], type="R", center=ER, limits=c(LCLRa,UCLRa))}
\end{verbatim}

Deze grafieken laten zien dat het gemiddelde van de lengte van de spuitjes duidelijk out-of-control is, het gemidddelde is veel groter dan gewenst.\\
--> vreemd is dat de variatie wel kleiner is dan `gewenst'. HOE KOMT DIT????\\
Een technicus werd er bij geroepen om de machine weer goed af te stellen, deze heeft dit twee keer gedaan maar beide keren resulteerde dit bij de volgende meting in gemiddeldes boven de upper control limits. Deze metingen zijn zichtbaar in de $31^e$ en de $32^e$ rij van de dataset. Dit zijn ook de laatste twee punten in Figuur \ref{fig:control2}.\\
Nadat de technicus de machine een derde keer opnieuw had afgesteld bleef de volgende meting binnen de control limits in de $\overline X$ chart en in de$R$ chart. Hierna werd weer doorgegaan met metingen nemen elke 15 minuten.\\


\section{Production Run and the Persistent Statistical Quality Specialist}
Nadat weer 15 metingen waren gedaan werden deze weer geplot met de control limits van de eerste 15 data punten. Dit is zichtbaar in Figuur  \ref{fig:control3}.
\begin{figure*}[hbtp]
\includegraphics[width=7cm]{paper-xbar3.jpg} \includegraphics[width=7cm]{paper-R3.jpg}\\
\caption{Control charts}\label{fig:control3}
\end{figure*}
\begin{verbatim}{
paperqccX3 <- qcc(paper2[33:47,2:6], type="xbar", center=mu, limits=c(LCLXa,UCLXa))
paperqccR3 <- qcc(paper2[33:47,2:6], type="R", center=ER, limits=c(LCLRa,UCLRa))}
\end{verbatim}
Over deze charts heeft een statisticus de volgende conclusies getrokken: de originele 15 punten die gebruikt werden om de $\overline X$ en de $R$ charts te defini\"eren waren niet in statistical control. De laatste 15 punten lieten ook een proces zien dat niet statistisch in-control was.
Maar wat hem ook opviel is dat nu wel de variatie gereduceerd was (in de laatste chart), en dit is alleen maar positief.\\



\section{Wat was er aan de hand?}

Het bleek dat de technicus de eerste twee keer dat hij de machine opnieuw moest instellen hij de instellingen voor de lengte naar beneden schroefde. Maar dit kon hij niet goed vastzetten, daardoor liet dit steeds los. De derde keer maakte hij hem hoger vast, dan werden de spuitjes wel te lang maar dan zou het wel blijven zitten.\\
Het bleek dat bij het `te laag' vastzetten, de grip zo slecht was dat het mee ging trillen met de machine wat er voor zorgde dat het gemiddelde omhoog ging. Wanneer de machine in een hogere stand werd vastgezet, kon alles wel goed vast blijven zitten. De gemiddelde lengte werd wel iets langer, maar de variatie was in zoverre gereduceerd dat er toch acceptabele spuitjes uitkwamen.\\
In een capability study van de laatste 15 punten is het duidelijk dat de  gemiddelde lengte wel iets meer werd, maar de variantie minder. Dit is zichtbaar in Figuur \ref{fig:finalcap}.\\
\\
Uiteindelijk werd er voor gekozen om dit deel van de machine te vervangen, de kosten hiervoor liepen goed op doordat dit onderdeel erg duur was en het proces moest worden stilgelegd, maar dit hebben ze naderhand allemaal terugverdiend doordat het proces hierna veel beter liep.
\begin{figure*}[hbtp]
\begin{center}
\includegraphics[width=7cm]{paper-PCAlaatste15.jpg}
\caption{Final capability analysis}\label{fig:finalcap}
\end{center}
\end{figure*}

\section{Conclusie}
Deze studie laat zien dat in praktijksituaties er vaak veel factoren mee spelen. \\
\\
\noindent
Ook blijkt dat in de praktijk lang niet altijd op een goede manier met de data om gegaan wordt maar dat er alsnog goede conclusies uit de data gehaald kunnen worden. 
Het is belangrijk dat je een goede data-inspectie doet, dat je de aannames van je capabilitystudie test en of je data in-control zijn. \\
Het is namelijk erg logisch dat je SPC-procedure niet werkt omdat de Control Charts niet gebaseerd zijn op goede analyses. \\
\\
\noindent
Als er eerder en of beter op de data gereageerd was, was de fout in de machine eerder ontdekt. Al bij het plotten van de eerste 15 samples was de fout zichtbaar, en had de technicus erbij geroepen kunnen worden. Maar ook op deze manier werd de fout ontdekt, en daar gaat het uiteindelijk om.\\
\\
\noindent
Het was wel leuk op SPC een keer `in de praktijk' toegepast te zien worden. Het geeft een beeld van wat er allemaal bij komt kijken en wat er allemaal verkeerd kan gaan. Ook geeft het een mooi voorbeeld van hoe iemand met `gevoel' voor het proces (de technicus) een statistisch erg gekke aanpassing aan de machine doet (de instelling van de spuitjeslengte omhoog doen terwijl het gemiddelde te hoog ligt).  Het is goed om met allerlei aspecten rekening te houden en goed in de gaten te houden wie er wanneer en waarom allemaal aan de machine sleutelt. Hieruit bleek maar weer dat in de praktijk ook niet statistici met het proces omgaan die dus door hun niet-statistische blik op de dingen bepaalde verbanden niet leggen. 


\end{document}